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Abstract 

We present hierarchical Bayesian methodology to perform spatio-temporal change of support 
(COS) for survey data with Gaussian sampling errors. This methodology is motivated by the 
American Community Survey (ACS), which is an ongoing survey administered by the U.S. Cen¬ 
sus Bureau that provides timely information on several key demographic variables. The ACS has 
published 1-year, 3-year, and 5-year period-estimates, and margins of errors, for demographic 
and socio-economic variables recorded over predefined geographies. The spatio-temporal COS 
methodology considered here provides data users with a way to estimate ACS variables on cus¬ 
tomized geographies and time periods, while accounting for sampling errors. Additionally, 3-year 
ACS period estimates are to be discontinued, and this methodology can provide predictions of 
ACS variables for 3-year periods given the available period estimates. The methodology is based 
on a spatio-temporal mixed effects model with a low-dimensional spatio-temporal basis function 
representation, which provides multi-resolution estimates through basis function aggregation in 
space and time. This methodology includes a novel parameterization that uses a target dynami¬ 
cal process and recently proposed parsimonious Moran’s I propagator structures. Our approach is 
demonstrated through two applications using public-use ACS estimates, and is shown to produce 
good predictions on a holdout set of 3-year period estimates. 
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1 Introduction 


The American Community Survey (ACS) is an ongoing survey that releases data annually, provid¬ 
ing communities with the current information needed to plan investments and services. The ACS 
was designed to produce reliable annually updated estimates for topics only previously available 
once every decade from the decennial census long-form, such as detailed demographic housing 
and socioeconomic topics. The structure of the publicly released ACS data is unique, consisting of 
rolling multi-year estimates (MYEs). In August 2006, the U.S. Census Bureau released its first set 
of one-year estimates for areas with populations greater than 65,000. Subsequently, in December 
2008, the first set of three-year estimates were released for areas with populations greater than 
20,000. Completion of three types of releases culminated in 2010 with the release of the five-year 
estimates for all standard tabulation areas including census tracts and block groups. It was recently 


fiscal year ( 

US Census Bureau, 

20U 

5a). For additional details on the structure of ACS, see US 

Census Bureau ( 

2015b 

) and 

Torrieri ( 

2007 

)• 


The shift from the decennial census long-form data to using MYEs from the ACS offers unique 
challenges and opportunities for data practitioners. One of these is to transform the survey esti¬ 
mates from one geography and/or time period to another in a manner that can account for sampling 
uncertainty (i.e., allowing for custom geographies and/or time periods). In general, procedures that 
allow one to perform statistical inference on a spatio-temporal support that differs from the support 
of the data (either in space, time, or space and time) is referred to as change of support (COS). In 
this setting, we let the survey data’s spatio-temporal support (e.g., census tracts and a 5-year time 
period) be denoted as the source support, and let the support of practical interest be called the target 
support (which might also be a standard census geography). In particular, we are interested in the 
problem where either the spatial and/or temporal support differs between the source and target. Al¬ 
lowing ACS data-users to obtain estimates and conduct inference on user-defined spatio-temporal 
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supports has recently been identified as an importan t problem by a National Academy of Sciences 


(NAS) panel ([National Academy of Sciences . 


2015b . 


Clearly, the introduction and removal of 3-year period estimates has produced a need for spatio- 
temporal COS among the ACS data-user community. Although there is interest in this particular 


COS problem, there are broader implications of 
statisticians and ACS data-users. For example, 


spatio-te mpora l 


COS that are of interest to federal 


McElroy (2009) proposed methodology for com¬ 


paring trends between counties that do not have compatible MYEs. A method of spatio-temporal 
COS would allow one to readily compare (time) trends across different geographies, with poten¬ 


tially increased precision. Another example can be found in 


Bradlcv et al 


(12014bh . where New 


York City’s Department of City Planning is interested in producing estimates of socio-economic 
variables on community districts (a geography not made available by ACS) 


Bazuin and Fraser, 

2013; 

Siordia, 

2013c 

b.a 


Beaghen et al 


2012 . 


2011 


a, among others) and/or time (e.g., see 


furthermore, under- 

Siordia et al.. 

2012; 

see 

McElroy, 

2009; 


, among others) is a consistent theme in the ACS literature. These re¬ 


occurring problems suggest that spatio-temporal COS would greatly enhance the utility of ACS 
period estimates. Thus, our primary goal is to introduce a method to perform spatio-temporal COS 
for survey data, with sampling errors that are reasonably modeled as Gaussian random variables. 
In the geographical sciences, a standard approach to the spatial-only COS problem is known 


of the target areal units (e.g.. 


Killough, 


1991 


Flo were 

ew anc 

Green. 

1989; 

Flowerdew et al.. 

1991 


Flowerdew and Green! 


1992. 


; Rogers and 


19941 . Although such procedures are easy to imple¬ 


ment, measures of uncertainty are not readily available, which limits their usefulness in situations 
with substantial measurement/sampling uncertainty. 

For reviews, 
). There are 


Spatial statistics is an avenue for spatial COS th 

at takes into account uncertainty. 

see 

Gotway and Young 

(2002 

), 

Cressie and Wikle 

(2011 

), and 

Baneriee et al. 

(2015 


two general methodological approaches currently used for spatial COS. The first method involves 
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defining a spatial process at the point level, and then integrating the process to the desired target 


support. For an example of this “bottom-up” strategy see 


Wikle and Berliner (|2005l) where they 


consider a hierarchical Bayesian approach for Gaussian data that includes simple areal interpola 


tion as a special case; also see 


Bradley et al. 


(I2014bh who consider a hierarchical spatial Poisson 


COS methodology. The second general modeling framework for spatial COS starts by defining 


parameters for each set in a partitioning between t 


this “top-down” strategy see 


Mugglin and Car lin 


le sou rce a nd target suppo r t. For examples of 


199811 and 


Mugglin et al. 


(11998b . where they 


use a hierarchical Bayesian approach based on a Poisson data model. In general, both of these 


approaches give similar results; see 


Trevisani and Gelfand (12013b for an explicit comparison. 


We choose to specify the latent process on a point-level spatial support and integrate to any 
desired target support (i.e., the “bottom-up” method for spatial COS). This allows one to avoid 
computationally expensive Bayesian simulation every time a new target support of interest is in¬ 


troduced (e.g., see the discussion in 


Bradley et al.. 


2015b 


c). This flexibility is especially valuable 


for our application, since ACS users may continually propose new target supports. 

Although the literature for spatial COS is mature, very little statistical work has been done 
on spatio-temporal COS. We focus on such a methodology for data that are reasonably mod¬ 
eled as Gaussian. The methodology explicitly accounts for sampling uncertainty in the survey 
estimates, as well as differing geographies and MYEs. We utilize a Bayesian hierarchical model¬ 
ing framework with a key assumption that the important spatio-temporal variability in the latent 
spatio-temporal variable of interest can be efficiently modeled in terms of a relatively low-rank 
spatio-temporal basis expansion, with associated random expansion coefficients. This modeling 
approach allows for the COS methodology to be applied to high-dimensional datasets. Although 
various basis expansion approaches have become the norm in spatial statisti cs ove r the last few 
years, both for spatial processes and spatio-temporal processes (e.g., see lCressie and Wikle . 


2011 


for an overview), these models have not typically considered joint spatio-temporal basis function 
expansions. The joint spatio-temporal basis functions are important in our methodology as they 
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allow for simple aggregation across spatial scales and time. 


In addition, we present a novel parameterization of the random effects dependence structure 
that makes use of a “target” dynamical spatio-temporal process and what we have termed the 


Moran’s I propagator (see 


Bra dle y et al. 


2015ah . This is important as it respects the fact that 


there is an underlying dynamical process, but recognizes that we only need be concerned with the 
implied marginal dependence suggested by that process. This parameterization is further charac¬ 
terized by its extreme parsimony, thus providing a fairly low-dimensional parameter space for our 
Bayesian hierarchical model. We illustrate this methodology by considering the MYEs for ACS 
median household income and predict the 3-year MYE for 2013 (which was held out of our analy¬ 
sis) to provide a demonstration that our methodology provides a principled approach for predicting 
the discontinued 3-year MYEs for ACS data. Furthermore, a demonstration of simultaneous spa¬ 
tial and temporal COS is given, where 1-, 3-, and 5-year period estimates and counties are used as 
the source support, and 1-, 2-, 3-, and 4- year period estimates and “American Indian area/Alaska 
native area/Hawaiian home lands” are defined as the target support. 

The remainder of this paper proceeds as follows. The spatio-temporal COS methodology is 
outlined in Section 2, with the ACS example presented in Section 3. We conclude in Section 4 
with a brief discussion. 


2 Methodology 

The data of interest rely crucially on the sample design as well as the direct estimates and their 
estimated variances. For the sake of generality we present our approach in terms of MYEs. Note, 
for ease of notation we refer to a 1-year period estimate as a MYE. Let z\ ( \a) represent the ACS 
data quantity of interest for the Ath time period (/:' = 1,3,5), the t -th year, and for spatial region 
A e D) A , where the set D) A is a collection of areal units that can be different depending on the 
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period and time. We then specify the following data model 

Z^\a)=Y^\a) + e^\a)-, A£D\% t = T l ,...,Tu, £ — 1,3,5, (1) 

where (A) is considered to be the “true” (but unknown) latent variable of interest at time period 
£, time t, and spatial region A, and the “survey error” £, ^ (A) is assumed to be an independent 
and identically distributed (i.i.d.) Gaussian random variable with mean zero and variance Of f [A). 
These variances are assumed known and computed from ACS estimates of the margin of error 
associated with (Z^(A)}. 


2.1 Process Model 


One of the challenges with spatio-temporal COS is to formulate a flexible spatio-temporal model 
for the latent process, Y^\a). It is intuitive, and customary in gcostati sties, to consid er COS 


Cressie and Wi kle. 


2011 . 


from the perspective of aggregating a point-scale latent process (e.g., see 
Section 4.1.3). That is, for any spatial region S CDs G W ! and period £, we can define the discrete¬ 
time, continuous-space, spatio-temporal aggregate as 


Y t [( \s)= £ ~f Y(s;k)ds, \S\ > 0, (2) 

k=Ti -i l 5 l J °s 

where \S\ — Id s ds, (Of are temporal weights, and Y (s ;k) is a continuous-space, discrete-time spatio- 
temporal process defined at spatial locations s G Ds and times k e {7} ..... 7)/ }. 

The aforementioned methodology relies on the specification of the “point-level” spatio-temporal 
process. In particular, consider the point-level spatial-temporal representation 


F(u;k) = 5(u) + Y, Yj( u ’ k ) Ij, ueD s ,keD t , (3) 

7=i 

where <5(u) is a large-scale spatial trend term, {t//)(u;k)}J =1 corresponds to a pre-specified set of 
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basis functions indexed in space (u) and time ( k ), and {rjy} are the random expansion coefficients 
associated with each basis function, assumed to be mean zero Gaussian and potentially dependent 
(see below). At this point, we leave the specification of the basis functions general, but present a 
specific example in Section 2.2. 

Upon substitution of © into ©, we have the following model for (A): 

y ” (A) = w\L 5{u)du+ m k J‘J^ l ' l ' liu '’ k),]idu ’ 

= P(A) + J £ Ij + UV), (4) 

1 k=t-£+ 1 j= 1 

where we have assumed co t = l/£, defined /i(A) = -^J A <5(u)du, y/j(A;k) = J A y/j(u;k)du, and 
we have assumed that, given a sufficiently large r and upon integration, the finite truncation of the 
basis expansion is a reasonable approximation. That is, we have assumed 

00 1 t l r 

L 7 L 777 / Vj(w,k)du r/j&O; A cR d ,t = T L ,... ,Tu,£ = 1,3,5. 

j=r+ 1 1 k=t-e +1 Ja 

To account for this truncation and potential model misspecification, we include the “fine scale” 
error term {^\a)}, which are assumed to be i.i.d. Gaussian random variables with mean zero 
and variance oj. Note that we are assuming this truncation error applies at the aggregate spatial 
scale. 

In practice, it is convenient to define a spatial support at which parameters are assumed to be 
constant at lower resolutions. Let Dg = {Bj : i = 1,... ,ng} be this pre-defined fine-level spatial 
support consisting of disjoint areal units. For example, in our application (Section 3), we let Dg 
consist of every US county as defined in 2013 (note that county definitions change from year-to- 
year). Now, let <5(u) = jli for any u £ W l such that u £ Bj, Bj is the z-th areal unit in Dg, /i, £ M 
is unknown, and i — \.... .ng. This implies that the spatial trend term is constant within each of 
the i = 1,... ,ng areal units in Dg (with the respective value /i ; ). If we partition A £ D y t ^ into its 
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potential overlap with all B e Dg then 


if if r , . S Ans, 

^) = W L J 5 ("M"= w E^y = 

' AnS; ' AnS; 


(5) 


Letting /z(A,/) = |Afll?;|/|A| and defining the n^-dimensional vectors h(A) = (h(A, 1),... ,/i(A,ns))' 
and ju B = (jui,..., !in B )', we have 


^( A ) — 777/ 5 (u)Ju = h(A)> fi ; A c M 2 . (6) 

Our final process model at the aggregate level of spatial and temporal support can then be 
written 

F/V) = h(A)' J u s + ^ ) (A)'r 7 +#(A), ( 7 ) 

where y/j^(A) = (y/^ (A),..., y/ff (A))' and y/^(A) = |Ea- =/ _^ + i Vj{A',k) are the spatio-temporal 
aggregated basis functions with {^/^(-)} denoting mean zero Gaussian i.i.d errors that are inde¬ 
pendent in time and space. The model in ([7]) is nonstationary, nonseparable, asymmetric, and 
provides a way to easily model areal units A C on different scales. To our knowledge, there 
has been no such model proposed that allows for realistic spatio-temporal covariances (i.e., nonsta¬ 
tionary, nonseparable, and asymmetric) that is flexible enough to model data defined on multiple 
spatial and temporal scales. 


2.1.1 Random Effects Parameterization 

The most critical component of the process model in © is the r-dimensional random effects vec¬ 
tor, rj ~ Gaw(0, K), where K is an r x r covariance matrix. Although one could specify a general 
Bayesian covariance prior (e.g., through an inverse-Wishart distribution, modified Cholesky de¬ 
composition, etc.), such a prior does not reflect the underlying spatio-temporal dynamical process 
occurring at fine resolutions. This motivates us to consider a parameterization for K that incorpo- 
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rates fine-scale (e.g., the B-scale) dynamics within the model stated in Section 2.1. Specifically, we 
define a “target process” that is dynamical on a single small-scale geography (the B-scale). Then, 
K is chosen so that the B-scale covariances of {(■)} in 0 are close to the analogous covari¬ 
ances of the target process, which in the Gaussian setting, incorporates dynamics within 
on the B-scale. This approach has important implications for modeling (temporal) multi-scale pro¬ 
cesses in that current parameterizations available for a single scale spatio-temporal process can be 
readily used to define a parameterization for multi-scale spatio-temporal processes. 

Let the “target-process” on the finest areal spatial scale be given by 


y/(B) = / 1 (B) + v,(B), t = t l ,...,Tu,b e d b , 


( 8 ) 


where /t(B) is defined in ©, and v, (B) is a Gaussian random variable with mean zero where 
V/ = (vAB) : B F Dr]' is the associated Gaussian random vector. As discussed in Cressie and 


Wikle (2 0111) . many realistic dynamical spatio-temporal processes can be represented by fairly 


simple first-order Markov models. Thus, assume that v t follows a first-order vector autoregressive 
model 

v ? = Mv f +b ; , t — Ti,... : Tu, (9) 

where M is a r x r re al-va lued propagator matrix. We make use of what we have termed a “Moran’s 


I” (MI) propagator ( Bradley et al. 


2015ah . which only requires knowledge of the adjacency struc¬ 


ture of our areal spatial domain. In addition, the MI propagator is helpful in that it can accommo¬ 
date parsimonious time-varying behavior through the use of spatio-temporal covariates. Thus, we 
set M equal to the MI propagator matrix (see the Appendix for the definition of MI propagator). 

Given M, let b ? be an 775 -dimensional Gaussian random vector with mean zero and precision 
matrix E^ 1 = (1 /<j£) (I — A), where > 0 is unknown, and A is the adjacency matrix formed by 
the edges implied by Z)g. Then, assuming stationarity, the VAR(l) structure suggests the following 
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(e.g., see lCressie and Wikld.l201ll Chap. 6): 


cov(Y,* +t ,Yn=4 t , ) =M’-y.> 


40) 


where 

vec(E^) = [I-M®M]- 1 vec(E & ). 
Thus, given o^, M and E&, one can specify all of the elements of 

Ey* = var{(Yj^,..., Y^/}, 

where Y* = ( Y*(A ) : B e D B )' for t — T L ,..., T/j. Finally, we let 


K = arg min{||E y *-'f'C'f' / || F }, 

c 


( 10 ) 


where > 0 is unknown, for generic square matrix G we have 11G11= trace (G / G) is the 
Frobenius norm, the space of C in (flQl) is restricted to be positive semi-definite matrices, and 
the £ n[ 1 * x r matrix 

'p=(4) ) (-4)',---,v'4w')'. (11) 


Notice that the period i is set equal to 1 in (fill) , as the target process is defined on t 


le finest spatio 


temporal reso 

ution. 

"he solution to (fTOl) is known and easv to compute fe.a.. see 

Higham, 

1988 

Bradley et al., 

O 

<N 

2015a.c 

; Burden et al., 

2015 

). 




The r x r covariance matrix K is extremely general and can easily be adapted to account for 
other spatio-temporal covariance structures in the multi-scale setting. That is, if one has a space- 
time covariance matrix (say E*) that is not readily interpretable on multiple spatial or temporal 
scales, then replace Ey* with E* in (fTOl) . This would produce a process {t/^(-)} in ©, that has a 
space-time covariance function that is close (in Frobenius norm) to E* on the scale in which E* is 
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defined. 


In summary, given the spatio-temporal aggregated basis functions, a MI propagator matrix, 
fine-scale adjacency matrix, and a single variance parameter, we are able to obtain in (fTOl) a positive 
definite marginal covariance matrix for our spatio-temporal random effects, tj , that respect a target 
dynamical process at the fine (aggregated) spatial scale. It is important to note that we do not 
predict v t nor Y*, they are just an intermediate step toward obtaining the target marginal covariance 


in order to apply the 


Highaml ( 1988 ) best positive approximate procedure. 


2.2 Parameter and Basis Function Specification 

The model specified above is extremely parsimonious for a spatio-temporal model. In particular, 
we need only specify oj, oi and p B . In the application presented in Section 3, we assign both 
a| and a/G(l, 1) (inverse-gamma) prior and the elements of /i are assumed a priori to be i.i.d. 
N(0 ,ctj), with <7^ assigned an 7G( 1,1) prior. 

The primary reason that our model can easily adapt to different levels of spatio-temporal sup¬ 
port is because we model the spatio-temporal variability through a basis expansion. In particular, 
we utilize spatio-temporal basis functions and common random effects. In this regard, it is im¬ 
portant to specify a basis set that is flexible and can easily be integrated/ag gregate d in space/time. 


There are many such choices in the literature (Wikle, 
sider local bisquare basis functions given by 


20101 : 


Bradley et al. 


2014c). Here, we con- 


Yj( = 


{i 

0, 


u-c, 


w s) ~(\t-gt\/wtj } , if ||ll — Cy|| <W s ,\t~g t \ < W t 

otherwise; u £ D s , 


( 12 ) 


with j = 1..... r spatial knot points Cy, m, equally spaced temporal knot points, g t , where w s and 
w t are maximal spatial and temporal supp orts. The placement of knots is chosen using a space 


filling design (N vchka and Saltzman , 


1998 ). Additionally, direct Monte Carlo sampling is used to 


approximate tg 7 (A;k). Specifically, h points {s g : q — 1..... fi} C A are randomly selected using 
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a uniform distribution on A. Then, t \fj{A\k) is approximated with (l//i)Lg=i |A| Yj(s q ;k), where 
Yj(s q ;k) is computed according to (fl2l) . In the ACS analysis considered in Section 3, we specify 
the components of the basis functions (i.e., w s , w t , and r in (fl2l) l. 


3 ACS Median Household Income Example 


The US Census Bureau has amassed a large number of ACS period estimates, and has become an 
extremely rich data source. In fact, a recent Google Scholar search (on ) of “American Community 
Survey” resulted in 3,690,000 entries. Thus, providing ACS users the flexibility to undergo space- 
time COS is likely to have a large impact. To demonstrate this, in Section 3.1 we perform space- 
time COS to show that one can provide estimates of the 3-year period median household income 
using the available 1-, 3-, and 5-year period estimates. This is e speci ally notable considering 


that ACS has decided to discontinue the 3-year period estimates (lUS Census Bureau . 


2015ah . 


Additionally, we provide an example of simultaneous spatial and temporal COS in Section 3.2. 


3.1 Estimating 3-Year Period Estimates of Median Household Income 

We consider using ACS period estimates of median household income defined at the county level. 
Modeling at the county level for this illustration is partly based on the broad recognizability of 
this level of geography to both those within and outside of federal statistics. Also, the spatio- 
temporal coverage of ACS period estimates defined on counties is comparatively better than other 
more spatially sparse geographies; for example, census tracts have no 1-year period estimates and 
relatively few 3-year period estimates of median household income currently available. 

We consider 1-year ACS estimates of median household income for 2006; 1-, 3-year period 
estimates for years 2007 and 2008; 1-, 3-, and 5-year estimates for 2009 through 2012; and 1- and 
5-year estimates for 2013. We exclude the 3-year period estimates for 2013 in the model fitting 
and estimation stage in order to keep them as a hold-out sample. To summarize, we define the data 
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for the periods 


Z ( /\a) :AeD^l (t,£) = (2006,1),..., (2013,1), (2007,3),..., (2012,3)(2009,5),..., (2013,5). 

(13) 

In total, there are 19 time periods considered here, which gives a large dataset of 32,836 ACS 
period estimates. As an example, Figure 1 shows the 2013 1-, 3-, and 5-year estimates of median 
household income along with the associated standard deviations (based on the ACS margin-of- 
error estimates). 

To compute the predictor 



(A) 


r/V)|{z, ( V)} 


AeD 


(0 

t,A’> 


t = Tl,...,Tu, £—1,3,5, 


where the posterior expectation is taken using the model introduced in Section 2, we need to specify 
the spatio-temporal basis functions. This includes defining w s (radius of the spatial component of 
the basis functions), w t (radius of the temporal component of the basis function), the number of 
knots, and the knot locations. 

We set w s to be 1.1 times the smallest distance between two different knots, g t is set equal to the 
mi mid-points associated with the different ACS period estimates of median household income, and 
w t — 1.1. The value for w s and w t were chosen to minimize hold-out error (i.e., squared distances 
between predictions and (■)}). Similarly, we considered r — 50,100,150,200,250 and 300 
and the value of r that minimizes hold-out error is r — 250. Thus, in the following application r — 
250 and t = 2005.5,2006,... ,2012.5 (i.e., m t = 19), so that t $\a) is a rxm t — 4750-dimensional 
vector. 

(3) 

In the fourth row of Figure 1 we display the predicted median household income ^ 2013 (')• 

Visually, these predictions perform quite well when comparing to the hold-out ACS sample (i.e., 

(3) 

{Z 2013 (■)}• Also, notice that the posterior standard deviation is considerably smaller than the ACS 
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standard deviation (the model based errors are close to 20 times smaller than the ACS errors). 
Figure 2 provides a second look at the hold-out versus the model-based predictions, and again, our 
predictions appear to be performing quite well. 

As a diagnostic measure, consider the ratio 


m = 


zS!,(a) 


2013 


?( 3 ) 

'2013 


(A) 


• A G D ^ 

, /i t ^2013,A' 


(14) 


Here, if R is close to 1 then the model-based estimate is similar to the hold-out ACS estimate, if 
R is less than 1 then the model-based estimate is larger than the hold-out ACS estimate, and if 
R is greater than 1 then the model-based estimate is smaller than the hold-out ACS estimate. In 
Figure 3(a), we plot a histogram of the values in the set {/?(A)}. The histogram is slightly skewed 
right, indicating that there is more of a tendency to under-predict than over-predict; however, a 
majority of the mass of the histogram is located at 1 (similar behavior can be seen in Figure 2). 
Thus, we see that we are consistently reproducing similar estimates to the ACS hold-out 3-year 

period estimates. Finally, to further corroborate the results using the ratio, in Figures 3(b) and (c), 

(3) 

we plot the histogram of R(A) associated with the hold-out experiment that sets aside Z 2013 (-) and 

(3) 

^2012 (')’ and obtain similar results. 


3.2 Example of Simultaneous Spatial and Temporal COS 


In 1975 Congress appended Section 203 to the Voting Rights Act, which provides voting resources 


for US citizens that are not proficient in English. Recently, 


Joyce et al. 


(2014) described a pre¬ 


cise approach to classify regions that satisfy the jurisdiction rule laid out by Section 203. Their 
results indicated that many American Indian areas/Alaska native areas/Hawaiian home lands, met 
the jurisdiction rule of Section 203. Thus, considering the need for language assistance among 
US citizens that reside in these areas, it would be worthwhile to determine whether additional 
assistance is required based on (low) income status. Note that it has been shown that those at 
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lower income_ levels tend to face certain obstacles (e.g., transportation needs), making it difficult 


to vote (see 


Gelman. 


2008 


, for more discussion). We use this important example to demonstrate 


simultaneous spatial and temporal COS. 

Using replicates from the posterior distribution in Section 3.1 we produce model-based 3-year 
period estimates of the median income of individuals that reside in American Indian areas/Alaska 
native areas/Hawaiian home lands. In Figure 3(d), we plot a histogram of the ratio in (fl4l) us¬ 
ing ACS 3-year period estimates of median income in American Indian areas/Alaska native ar¬ 
eas/Hawaiian home lands. As in Figure 3(a,b,c), the ratios are consistently close 1 indicating 
strong out-of-sample performance of our proposed method. 

In Figure 4, we provide time series plots for 3 regions in our target support: the Navajo Nation 
reservation and off-reservation trust land, the Uintah and Ouray reservation and off-reservation 
trust land, and the Wind River reservation and off-reservation trust land. This is done to highlight 
the ability of our method to compare across geographic regions over time (that are incompatible 


when using ACS c 
(e.g., see McElrov 


ata alone), which is a problem of interest among the federal statistics community 


2009) . We choose these specific regions because they display the most notable 


patterns in median income, and in general, the time series plots associated with each region display 
different patterns. In Figure 4 (a,c,e,g), we see that the Navajo Nation reservation has low median 
income, and that both the Uintah and Ouray reservation and the Wind River reservation have 
median incomes that increase and then decrease over time. Then in Figure 4(b,d,f,h), we see that 
we obtain precise estimates using our approach (i.e., the posterior standard deviations are small 
relative to the scale of the data). Furthermore, the standard deviations appear to be larger for earlier 
years than for later years, which conforms to intuition since more ACS estimates are available as 
time goes on. Additionally, as we increase the period we see that the posterior standard deviations 
are smaller. 
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4 Discussion 


In its relatively short period of existence, the ACS has proven to be a valuable resource for data 
users across academia, government, and industry. In many cases, data users prefer to have esti¬ 
mates at different geographies and/or time periods than are provided by the U.S. Census Bureau. 
We present a novel approach to spatio-temporal COS that allows data users to take published ACS 
MYEs and provide predictions at arbitrary spatio-temporal levels of support. Our approach is novel 
in that it is based on a spatio-temporal mixed-effects model that relies on a spatio-temporal basis 
expansion, wherein the basis functions are easily aggregated and the associated random effects are 
independent of scale. In addition, we provide a novel and extremely parsimonious dynamic model 
motivation for the marginal covariance matrix of these random effects. This representation allows 
implementation on extremely large datasets. We illustrate the effectiveness of our methodology on 
a holdout sample of 3-year ACS MYEs for 2013, analogous to a potential real-world implemen¬ 
tation that will arise when 3-year MYEs are discontinued in 2016. Additionally, we demonstrate 
simultaneous spatial and temporal COS from ACS county-level period estimates to 1-, 2-, 3-, and 
4-year model based estimates of median income defined on American Indian areas/Alaska native 
areas/Hawaiian home lands. 

The approach proposed here has clear utility beyond the ACS application presented in Section 
3, and even outside of the small area estimation context of federal statistics. Specifically, with a dif¬ 
ferent choice of basis functions and possibly a different dynamic structure, this methodology can be 
applied in areas as diverse as environmental science, and ecological modeling among others. Nev¬ 
ertheless, the potential impact to data-users and policy makers interested in ACS custom-designed 
tabulations is unparalleled. 

The approach we propose provides several opportunities for methodological extension. In par¬ 
ticular, count data is prevalent within the ACS and other federal surveys; thus, there is scope for 
extending spatio-temporal COS for non-Gaussian data, as was done in the spatial only case in 
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Bradley et al. 


(1201461) . Additionally, to date there is relatively little (by comparison) literature in 


the context of multivariate spatial and spatio-temporal COS. The methods proposed here and in 


Bradley et ali (I2015ah provide the initial building blocks for developing a comprehensive frame¬ 


work and is a subject of future research. 


Appendix: Review of Moran’s I (MI) Propagator 


We_provide a review of the MI propagator matrix, which was recently introduced in 


Bradley et al. 


(l2015ah . Confounding in mixed effects models is the core motivator for the MI propagator matrix, 


but it has additional benefits related to model parsimony and the ability to accommodate covariate- 
based non-autonomous propagators. By confounding, we mean that the columns of the design 
matrix are linearly dependent with the columns of the coefficients of random effects. 

To see this, rewrite ([8]) in vector form as, 


Yt = H,u + v,, 


(A.l) 


where H = (h(A)': B e D B )' and t — T Ll . ,.,Ty. Then, substitute © into (IA.il) to obtain 

Y* = H/i + M ; jU f _ j + b,; t = 7^,...,7)/. (A.2) 

Depending on our choice for {M f } there might be issues with confounding between V t -\ and the 
2n-dimensional random vector Q t = (/T, bj)'; t — T L ,..., Ty. Upon rewriting (IA.2I) . we get 

y;=BC, + M,v,_ i; (A.3) 

where the n x 2 n matrix B = (H, I). The strategy, is to set the columns of the propagator matrix M f 
equal to columns in the orthogonal complement of the column space of B. This ensures that the 


16 













columns of B are linearly independent of the columns of M f . Now, using the spectral representation 
of (I - B(B / B) 1 B')W(I - B^'B^B') = we set the r x r real matrix M, equal to the 


first r columns of <t>g for each t, which is denoted with Mg. Here, any rxr re al-valued matrix W 


can be used, and in Section 3 we set W = I as is done in 


Bradley et al. 


(2015a). 
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(e) 2013 1-year ACS Estimates 




(f) 2013 1-year ACS Estimates of Std.Dev 
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Figure 1: The first three rows present 2013 ACS period estimates of median household income 
and their respective standard deviations. The period is indicated in the title. White locations 
are missing. The last row presents our model based 3-year period estimate of median household 
income over all counties along with the associated posterior standard deviations. 
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Figure 2: Predictions based on all data except {Z}g 13 (A)} is given by the black lines. The estimates 
in the set {Z^, 1 3 (A )} are given by the red lines. 
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(a) Histogram of the Ratio of ACS to Model Based Estimates 
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(c) Histogram of the Ratio of 2012 ACS Hold-Out to Model Based Estimates 
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(b) Histogram of the Ratio of 2013 ACS Hold-Out to Model Based Estimates 
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(d) Histogram of the Ratio of ACS to Model Based Estimates 
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(3) 

Figure 3: The ratio in (5). Panel (a) gives a histogram of the ratio between {Z 20 | 3 (A)} and the 

(3) (3) 

predictions that hold out {Z 20I3 (A)}. Panel (b) gives a histogram of the ratio between {Z 20P (A)} 

(3) (3) 

and the predictions that hold out both {Z 2013 (A)} and (Z 20p (A)}. Panel (c) gives a histogram of 

(3) (3) (3) 

the ratio between {Z 20p (A)} and the predictions that hold out both {Z 20L3 (A)} and (Z 2012 (A)}. 

Panel (d) gives a histogram of the ratio between the 2013 3-year period ACS estimates of median 
income over American Indian area/Alaska native area/Hawaiian home lands, and the corresponding 
model-based estimates (see Section 3.2). 
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Figure 4: The first column gives time series plots for 1-, 2-, 3-, and 4-year model based period esti¬ 
mates of median income for the Navajo Nation reservation, the Uintah and Ouray reservation, and 
the Wind River reservation, respectively. The right column displays the corresponding posterior 
standard deviation associated with the estimates given in the first column. The legend indicating 
the American Indian reservation is given in Panel (h). Notice that the each row has a different 
range of years indicated on the x-axis, and the y-axis differs in Panels (b), (d), (f), and (h). 
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